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We integrate the MHD ideal equations of a slender flux tube to simulate the internal plasma 
' dynamics of coronal post-flare loops. We study the onset and evolution of the internal plasma in- 

stability to compare with observations and to gain insight into physical processes and characteristic 
parameters associated with flaring events. The numerical approach uses a finite-volume Ifarten-Yee 
TVD scheme to integrate the fD| MHD equations specially designed to capture supersonic flow 
\ discontinuities. We could reproduce the observational sliding down and upwardly propagating of 

. brightening features along magnetic threads of an event occurred on October fst, 200i. We show 

that high-speed downflow perturbations, usually interpreted as slow magnetoacoustic waves, could 
be better interpreted as slow magnetoacoustic shock waves. This result was obtained considering 
adiabaticity in the energy balance equation. However, a time-dependent forcing from the basis is 
O ' needed to reproduce the reiteration of the event which resembles observational patterns -commonly 

known as quasi-periodic pulsations (QPPs)- which are related with large scale characteristic longi- 
C/3 , tudes of coherence. This result reinforces the interpretation that the QPPs are a response to the 

pulsational flaring activity. 
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Oscillatory processes of several minutes periods in stellar coronae are associated with large scal e quas i-periodic 
. pulsations (QPPs) in flares. They were identified in X-ray stellar coronal flares (Mitra-Kraev et al. l2005f ) and were 
• ' later reported for the S un by FouUon et al. 2005j |. Shorter period QPPs in the radio band are known since the 
] early 70's f Aschwanden |2063 ). Also, systematic intensity perturbations in post-flare loops could indicate that they 



are the result or evaporation/ condensation, or rising/ lallme plasma now cycles, caused by the enicient heating ot the 
flaring plasma from the chromosphere (De Groof et al. l2004f ). This suggests that energy released by flaring pulsations 
can induce oscillatory loop trapped processes. Many clear quasi-periodic phenomenon are believed to be associated 
. ^ with wave responses after chromospheric evaporation has taken place as the main initial matter inflow source for 
^> . flare loops (Nakariakov et al. ,2004 ). Recently, Sigalotti et al. j2009i] simulated the propagation and damping of 
^ ' localized impulsive perturbations in coronal loops considering different dissipative mechanisms in the linear limit. As 
5^ ' the authors state there is still no consensus over the actual mechanisms from which the waves originate. Moreover, 
the importance of different damping mechanisms is not fully understand due to the fact that different theoretical 
scenarios can describe similar observational frames. In this paper we confine to the ideal MHD approach -allowing 
the consideration of nonlinear waves and shocks which give rise to initially non-dissipative damping mechanisms- to 
evaluate its capability to adjust observational scenarios. 

In Borgazzi and Costa j2005j we could describe an observational framework of catastrop hic eva cuation or high-s peed 



descending flow picture resembling other authors results e.g. Kjeldseth-Moe and Brekke jl998j |: Schrijver [2001 1. An 



upwardly propagating pattern was also registered associated to other parts o f the loop system. Also, the frequent 
recurrence of longitudinal compressible oscillations (e.g. in Borgazzi and Costa '2005], a quasi-periodic iteration of the 
event occurred 3 times in 1 h) were proposed to be ass ociated with cither large wavelength slow magnetoacoustic waves 
or limit-cycle solutions. In Costa and Gonzalez 2006j | we founded that limit cycle solutions associated with flow-based 



models are convergent with large wavelength wave-based models. Magnetoacoustic modes having wavelength of the 
order of the loop structure describe a basic oscillation between the kinetic energy of the parallel plasma flows and the 
internal energy i.e. the thermal perturbation. The structure of the obtained mode is such that in half of the period 
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the disturbance is always positive and the plasma emerges from the chromosphere while in the other half the plasma 
is negative and descending from the corona to the chromosphere. 

In observational scenarios the dynamics is complex making it difficult to isolate the erratic flaring from the oscilla- 
tions. Thus, due to the simplifications that are assumed for simulating physical quantities, the theoretical descriptions 
generally enables us to distinguish the behaviour of these components. Due to the hight Reynolds number of the corona 
it could be expected that a perturbation could give rise to a slow magnetoacoustic wave or a rising/falling mechanism 
that iterates before it decays in a non-dissipative way. 

Our intention is to provide a more accurate description of flaring multi-loop systems which exhibit a combined 
behaviour, sliding down and upwardly propagating appearance, to analyze the plasma dynamics interior to the loops 
and its relation with the flaring event, i.e. whether the oscillatory patterns -which could be associated with MHD waves 
and/or limit-cycle flow solutions in coronal structures- could be associated to internal or external driving mechanisms. 
As internal mechanism we refer to the development of non-dissipative second standing acoustic modes believed to be 
responsible for the induction of QPPs. These QPPs can be understood in ter ms of acoustic auto-oscillations, similar 
to the auto-oscillations in an electric-circuit generator (Tsiklauri et al. l2004[ ). The other alternative implies that the 
QPPs are the resulting features of the external forcing produced by the flaring event. 



II. NUMERICAL METHOD 



The basic MHD equations, the conservation of mass, conservation of momentum, energy balance, the magnetic 
induction equation together with a state equation, written in the so called conservative form result an hyperbolic- 
parabolic system of equations. The hyperbolic terms represent the convective effects and the parabolic ones the 
diffusive parts. These last terms are not considered in this paper so we work with an ideal MHD hyperbolic system 
of equations. We developed a ID i numerical technique that consists of an approximate Riemann solver that allows 
to evaluate the variables inside each cell by means of the variation of the flows through the contour of the cell. The 
software that implements this numerical technique, has proven to be robust solving very demanding benchmarks 
such as magnetohydrodynamic Riemann problems and Hartmann flows. The MHD Riemann problem solutions are 
comp ressible non-stea dy flow s; the Hartmann flows are incompressible and steady-state flows (Elaskar, et al. i200Ga; 
l2006b : Maglione et al. l2007t ). Thus, the technique is seen as a reliable computational tool for the description of MHD 
flows. 

Due to its accurate performance we apply the Harteen-Yee TVD technique (Total Variation Diminishing) which 
is specially designed to capture discontinuities when dealing with supersonic flows and to achieve a second order 
approach where the solution is smooth (Yce ct al. 1985). Thus, the one non-dimensional equations considering 
non-stationary flow and written in conservative form are expressed as: 

Ut + F(U), = (I) 

where U is the state vector of variables 

U = [p,pu^,pUy,pu^,Bx,By,B^,ef (2) 

and the hyperbolic fluxes are 

F = (pux,pul - BI+p + ^B'^,pu.j;Uy - B^By.pu^u^ - B^B^, , 



1 \T 

u^By ~ UyBx,UxB^ - UzB^, {e+p + -B'^ ju^ - B • uB^ j (3) 

p indicates de density; u^, Uy and Uz are the velocity components; B^, By and B^ represent the components of the 
magnetic field vector; p is the pressure and B^ = B^ + By + B'^. The total non-dimensional energy is 

^j + ^p{ul + ul + ul) + ^p{Bl + Bl + Bl) (4) 
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where 7 is the ratio of specific heats. 

In a quasi-linear form equation [1] results 



Ut Ae • U, = (5) 
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where represents the matrix associated with the Jacobian fluxes, letter 'c' indicates that the derivation is with 
respect to the conservative state variables. A simpler form of the Jacobian fluxes is obtain as a function of the 
primitive variables 



W = {p,Ux,Uy,u^,Bx,By,B^,py 



Equation [5] in primitive variables reads 



with the transformation rule 



Ap • W, = 



W, + A, • 



(6) 
(7) 
(8) 



The conservative form (eq. 5) is generally used for n umeric al purposes. They are appropriate to obtain accurate 
jump conditions at discontinuities and shocks (Leveque l2004[ ) and they ensure that mass, momentum and energy are 
conserved. However, the eigenvalue and eigenvector manipulation is simpler in the primitive form (eq. 7), thus we 
solve the MHD equations using the conservative form and the eigenvectors are evaluated by means of the primitive 
formulation. 

The resulting Jacobian flux has a null eigenvalue implying that is not possible to use a Ri eman n solver and thus 
an alternative Jacobian flux must be implemented. We used a technique developed by Powell jl995j | and we obtained 
the primitive eigenvalues which resulted 

Xe=Ux; Xa^UxCa] Xf^U^Cf] Xs ^ U^Cs] Xd ^ Ux 

where e stands for entropic, a for Alfven, and s and / for slow and fast respectively. The Alfven speed and the fast 
and slow speeds are respectively 
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We normalized the eigenvectors following the proposal of Zachary et al. 1994| to avoid vector degenerations. 
The explicit TVD scheme can be written as 



where the functions that determines the numerical fluxes are defined as 
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R is the matrix formed by the right eigenvectors of Ac. The dissipation function is expressed as 



(9) 



(10) 



(11) 



We a lso us ed a "limiter" function approapiated to optimize the one-dimensional flow (Elaskar et al. j200Cll | ; Maglione 
et al [2Q03]). 



9? = s5n(A™ Omax 0,min(a::;. | a™ . T^a™.) 



(12) 



with a{z) =1 z I if I z |> e and ^^(z'^ + e^) < e; ^j^Cgfli - 9"') if i 7^ and if a^^i = 0. 
The calculation process becomes simpler using the primitive variables to obtain 

is the left eigenvector of matrix Ap associated with wave m. 



(13) 
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FIG. 1: a) TRACE location of the brightening, measured along the loop, as a function of time. The spatial origin is chosen 
at the loop's Northern leg. b) Superposition of different events: locations of the brightenings as a function of time measured 
along a TRACE loop image: * for the first TRACE sub-interval; A for the second TRACE sub-interval and □ for the third 
TRACE sub-interval. To compare the two telescope results the MICA curve was superposed {<>) on the TRACE first case. 
The different lines of each time sub-interval are the asymptotic curve values. 



III. THE PROBLEM 

In Boreazzi and Costa [2005j we studied a post-flare event occurred October 1st, 2001 w hich was imagined by both, 
TRACE space telescope (Transition Region and C orona l Explorer, Handy et al. Il999l ) and MICA land telescope 
(Mirror Coronagraph for Argentina, Stenborg et al. Il999l ) . We described an observational framework of catastrophic 
evacuation or high-speed descending perturbations which resembles QPPs as it was repeated 3 times in 1 h. An 
upwardly propagating pattern was also registered associated to other parts of the loop system. Different hypothesis 
are given in literature to justify why a lternatively upward or downward evolving features are not seen (see the 
discussion in Borgazzi and Costa, |200.5| and references therein) . Figure [T^ shows the location of brightening as a 
function of time measured over a virtual axis that extends along a whole complex loop i.e. a compound system of 
neighboring coronal isolated threads which exhibit a coherent behaviour. The origin of the axis was taken at the 
loop's Northern foot and the distance was measured along the loop, considering that due to the location of the loop 
system and to the rough accuracy of the method projection effects could be discarded. The figure shows the sliding 
down of brightening features from the apex or bifurcation point towards both legs of the loop. Figure [T|d shows the 
superposition of events occurred in the same loop system at different sub-intervals of time. The coherence is due 
to the fact that the figure represents the spatial and time behaviour of several apparently isolated tubes. Moreover, 
as we found that there was a time correlation between relative maxima when comparing different structures of the 
same active region we suggested that the large longitude chromospheric coherence founded i.e.~ 300 Mm, must be 
associated with the forcing from the basis of the whole active region due to the intermittent flaring. However, this is 
not conclusive, recently Nakariakov et al. 2006j | developed a model that shows that QPPs observed in a flaring loop 



can be triggered by MHD oscillations in another loop situated nearby, not necessarily magnetically linked with the 
flaring one, and thus giving another possible explanation for the coupling of oscillations in an active region. 

Solar coronal conditions with large Reynolds numbers are well fitted by ideal MHD plasma models (i.e. infinite 
electrical conductivity cr 3> 1 leading to vanishing viscosity and ohmic dissipation) . We assumed that sources and sinks 
compensate each other, so the adiabatic energy equation holds, to investigate if the reiteration of features (Figure [Dd) 
can be associated with the internal dynamics of the plasma, i.e., the second spatial harmo nic of the acoustic mode 
which is believed to be responsible for QPP associated with fiaring events (Tsiklauri et al. l2004f ). These undamped 
long wavelength solutions can be also thought of as limit cycles. These are nearly longitudinal magnetoacoustic 
modes that describe a basic oscillation between parallel plasma kinetic energy and plasma internal energy with a 
characteristic length of L/4, L the length of the loop resulting in the balance between radiative losses and thermal 
conductive fiux (Costa and Gonzalez |2006| ). Otherwise, the reiteration of features could be associated to an external 
forcing from the chromospheric basis. This solution is generall y sugg ested by the chromospheric coherence observed in 
the dynamics of many coronal events (e.g. Borgazzi and Costal200i; Martinez et al. ,200a) . Of course, this description 
is not able to reproduce limit cycle models of the type proposed by Kuin and Martens |l982 | which are asso ciated 



with the generation and absorption of heat in open loop systems, neither the model proposed by MuUer et al. [2004 
via the assumption of different radiative loss functions, not taken into account here. 

Firstly our intention is to reproduce an individual observational scheme. And secondly we investigate if the reiter- 
ation of features could be explained by an internal triggered mechanisms or if it requires an external forcing one. In 
the frame of this model the role played by the magnetic field is mostly to guide the plasma fiow, so the simulations 
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FIG. 2: Scheme of the loop description. 




FIG. 3: Non-dimensional density -with respect to the background value- as a function of the non-dimensional loop length for 
different time steps, a) symmetric case, b) asymmetric case. The reference density and longitude values are po = W^^cm~^ 
and Lo = 150Mm respectively. 

are performed varying the thermodynamical quantities of the loop in appropiate coronal and chromospheric ranges. 

IV. RESULTS AND DISCUSION 
A. Simulation of an individual sliding down event 

We represent the loop system as a straight slender flux tube divided into three sections {xi, X2, X3), xi-\-X2+X3 = 1 
to integrate the non-dimensional ID i MHD ideal equations. The non-dimensional component of the magnetic field 
along the loop is = 1; By — 0.01, is an initial transverse perturbation, note that B^ >> By; the reference 
magnetic field value used is Bo — 20 G. The non-dimensional variables are taken with respect to the second section 
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FIG. 4: Simulation of the individual a) symmetric case and b) the asymmetric case. The observational data are superimposed 
.Dot lines correspond to the asymptotic observational values in Figure [T^. 
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FIG. 5: Non-dimensional density as a function of the non-dimensional loop length (same reference values as in Figure |3} for: 
a) time step 1000 (3.84 min). Shock wave fronts xis = 0.2; X2s = 0.8, contact discontinuities xic — 0.3; a;2c = 0.7, expansion 
waves xie = 0.45; a;2e = 0.55; b) time step 2000 (7.82 min). Shock wave xis = 0.05; X2s = 0.95, the expansion waves have 
interacted with the contact discontinuities and lowered the density and pressure. 



values, e.g. vi — ui/u2. Several simulations for different sets of initial conditions were performed. The density values 
chosen are pi = (10^, 10^°, 10^"'^) cm~^, the reference density value used is po = 10^^cm~^; the temperature values are 
Ti = (10", 10^, 10*^, 2 10*^, 10'^) K with a reference value of To = (10^ K). The pressure is obtained from the ideal state 
equation. The reference longitude value is Lo = ISOAfm. For different combinations of the state values we obtain the 
correspondent evolutions. Figure [2] gives a scheme of the initial conditions. 

Symmetric case 

We start with a symmetric partition to study the sensitivity of the system to the initial conditions, i.e. Xi = 
AQ%Lo,X2 = 20%Lo,X3 = AO%Lo. Table I displays the six initial condition cases that adjusted the upper slope data 
of Figure [T^. They correspond to a temperature step initial condition of the unique non-dimensional case, i.e., all the 
cases in Table 1. We describe the spontaneous evolution of the system from this static initial condition; the initial 
speed value used is Vi = 0. Cases 1 and 111 can be associated with cold models and thus short loops. In the frame of this 
model, the fact that the case that adjust observations results from initial constant density reinforces the hypothesis 
of homogeneous density distribution of most loops. Figure [3^ shows the evolution of the density front for different 
time steps. Figure is obtained from Figure [3^ assuming that the maximum of the perturbed density is a tracer 
of the observational features. We also obtain a resulting transverse magnetic field component smoothly modulated 
by the density front (the amplitude of the p ertur bed field is 10% of the non-perturbed one). However, following the 
classification by Nakariakov and Verwichte 2005j as the mode does not produce any noticeable perturbation of the 



loop minor radius it will be considered as a longitudinal mode, i.e., the magnetic field plays the role of being the wave 
guide of a fundamentally hydrodynamic shock. 

A more detailed description can be obtained from Figure [S^-b and from the wave analysis provided by the numerical 
techniques used. The figure displays the time steps number 1000 (3.84 min or 0.16 To To is the acoustic transit time 
along half of the loop) and 2000 (7.82 min or 0.32 Tq) respectively. Time step corresponds to p(0) = 1. FigurelH^ shows 
two descendant shock wave fronts {xis = 0.2; X2s — 0.8), two descendant contact discontinuities {xic = 0.3; X2c = 0.7) 
and two descendant expansion waves {xie — 0.45; X2e = 0.55). As the shock travels, the temperature, the density 
and magnetic field values are increased making the energy rise to values that could allow the observational detection. 
Note that the observation of moving brightening requires that these features must be time and spatially localized and 
thus, to make the detection possible, the density and temperature must be abruptly diminished again to compound 
a brightening feature. This effect is provided by the contact discontinuities which is recognized because the pressure 
and the velocity of the fiow are not changed while the wave passes. The expansion waves are recognized because they 
diminish the density and pressure values. From the comparison of time step 1000 and time step 2000 (cases a and b of 
Figure [5]), we see that the expansion wave has a larger speed than the contact discontinuity, i.e. when the expansion 
wave reaches the contact discontinuity the density and the pressure of the whole wave system is diminished. From 
Figure [3] we also note that the expansion waves have initially an ascendant phase (compare the density distribution 
at the loop center of steps 200 and 400 of the figure); the waves collide at the center and then travel downwards. 
Figure [Sh shows for time step 2000 the position of the descendant shock fronts {xis = 0.05; a;2s = 0.95). Also, the 
expansion waves have interacted with the contact discontinuities and lowered the density and pressure leaving behind 
a coupled nonlinear system of waves, i.e., the features located at xi = 0.2; 2:2 = 0.8 in Figure [5Jd can no longer be 
interpreted as contact discontinuities. 

Asymmetric case 

To simulate the asymmetry in the location of the bifurcation point of Figure [T^-b we assume xi — 0.31, X2 — 
0.12, X3 = 0.57, which is in correspondence with the asymmetrical dynamical geometry of the loop system, see 
Figure [T^. The aim is to reproduce the downward branch slopes of Figure [IJi, taking into account that to make the 
signal detectable the maxima must be higher than the background intensity value, e.g. Figure Starting from the 
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non-dimensional symmetric case, the pressure and density of the segments Xi are changed to reproduce the slopes. The 
slopes result strongly sensitive to variations of the density, i.e., note from Figure [Sb that higher background densities 
(right line starting at a; = 0.3) imply lower wave speeds (compare the speed of the right and left perturbations). We 
obtained few cases that could reproduce the conditions. The new non-dimensional initial condition values of segments 
{xi,X2,X3) that adjusted the observations were p (1.1,1.1,3) pressure p — (0.02,1,0.025) and T — (0.02,1,0.01). 
Figure HJd is obtained in the same way than in the symmetric case, assuming that the maximum of the perturbed 
density is a tracer of the observational features. The observational curve is superimposed over the numerical one 
showing good accordance between them. 

In the frame of this first model, it is not possible to give account for the reiteration of the brightening features of 
Figure [1)3. As a consequence of the low density values produced by the pass of the expansion wave, the temperature 
and density acquire very low values in the back front which are not recovered in the ascending phase. In all cases, 
the evolution from initial values systematically show a deep vacuum at the back front, i.e., density and temperature 
diminish to almost 30% of the initial background values. This can also be seen from Figure [3l the center values 
are lowered below the background density. However, to investigate if the low density region could be narrowed we 
vary the Xi partition. In this case, we obtain that the slope of the curve remained the same as in Figure and the 
vacuum region is narrowed but could not be lowered more than 8% of the total length, e.g., for X2 = 2 the low density 
region covers 8% of the center of the loop. Moreover, the fact that Figure [TJa shows only descendant features and 
not ascendant ones could be associated with the instrumental resolution that provides a limit value to the detection 
of the density. This implies that the density associated with the brightening of all the observed descendant branches 
must be higher than the ascendant ones, and last ones presu mably not detectable; this fact could not be reproduced 
by the numerical simulation of this model. In Costa et al. 2009j | we show that the pass of this type of non-linear 



waves form a vacuum tail that can give account of the so-called tadpole phenomenon. 



B. Simulation of the external forcing from the bases 



The different branch speeds of an individual phenomenon were adjusted changing the initial conditions, as in the 
first model approach given in subsection (4.1), but it was not possible to reproduce its iteration. To reproduce the 
sliding down features, as in Figure [TJa we decided to investigate the second triggering mechanism proposed. This 
is the chromospheric forcing from the loop footpoints produced by a transient oscillation associated with individual 
flare burst. The possible forcing functions are limited by the fact that they must reproduce the speeds and cadences 
of Figure [T|d. The temporally and spatially localized impulse deposition, resembling a mass injection in the bases of 
the loop, was switched off before the start of the simulation. Thus, the forcing must be considered as a boundary 
condition imposed to an homogeneous loop, i.e., with initial non-dimensional density p = (1, 1, 1) and temperature 
T = (0.1,0.1,0.1). Figure shows the density variation obtained assuming the action of one individual flare burst 
modeled by a forcing function f = A sin^((jjt), with A — 3, and the frequency u) = 6.28, which corresponds to a 
deposition of energy of ~ lO'^erg cm^^ on a typical loop. The value uj = 6.28 represents a dimensional characteristic 
time of T = 9.3 min. The Xi partition used is the same as in the first symmetric case. 

All the descendant branches have larger values of the perturbed density than the ascendant ones. The difference 
varies between ^ 20% and 40% (for the non-dimensional x value varying between 0.2 and 0.8) and thus, if the model 
is accurate, i.e., only the descendant branches are detected, an intensity threshold H can be numerically obtained. In 
fact, in the EUV and X-ray band, the emission intensity is modulated by the density perturbations as H oc {po + Sp)^ 

Ha = {po + 5p)\ <H<{po + Sp)l - Ho (14) 

and the threshold can be estimated by the above relation, po is the background density, dp the perturbed density, A 
and D means ascendant and descendant respectively. Taking into account the numerical results, i.e., non-dimensional 
densities (po + ^p)a = 2.45 and {po + 5p)o = 1.7, we obtain a mean value oi < H >= 4.446 10~^^ cm~^ with a 
AT? = 3.113 10~^^ cm~^ and AiJ/ < H >= 0.7001. Figure[S)3 is obtained from FigurelH^ using the intensity threshold 
\/ < H > — 2.11 to resemble the observational features, i.e., with only descendant branches. The ascendant branches 
are more intense than the descendant ones, however, extended ascendant features {0 < x < 0.2 and 0.8 < x < 1) can be 
distinguished at the location of chromospheric feet in Figure |61d. These features are correlated with the accumulation 
of observational points at the chromospheric footpoints in Figure [TJ) and could be indicating the wave front bouncing. 

The dimensional speeds of the numerical descendant branches are: vi — 106 km sec~^, V2 = 87 km sec~^, and 
V3 — 86 km sec^^, respectively. To compare with observation we take the mean value of the speeds of the northern 
and southern descendant branches of each event of Figure[TJ). The mean observational speeds are: ui = 129 km sec~^, 
U2 — 82 km sec~^ , and 1*3 = 65 km sec~^ , respectively. The time separations between branches for the numerical case 
are: Ti,2 ~ 24 mm, T2.S, 28 mm, and T^^^ 28 mm and the correspondent observational quantities (Figure [Dd) 
are: Ti,2 12 mm, T2.3 ~ 26 mm, and T^^^ ~ 22 min. 
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To compare the two model approaches, i.e., spontaneous and forced, we display the density profiles obtained for 
different fix temporal values of FigurelHl t = 2.23, t = 2.57 and t — 3.30. The resulting features are shown in Figure[71 
The density profile at t — 2.23 is shown in Figure [7^. We see two descending shock fronts at the center of the figure, 
i.e., xis — 0.4 and X2s = 0.6, together with ascending waves corresponding to the rebound of the first shock. This 
last features have their maximum at xi 0.1 and X2 ^ 0.9 respectively and both are below the non-dimensional 
density threshold corresponding to \/Hnd ~ Ph ^ 2.11. The density profile at t = 2.57 is shown in Figure [7)3. The 
two descending shocks are at xis ~ 0.2 and X2s ~ 0.8. At the center the nonlinear interaction of waves lowers the 
density. Only descendant shock features are detected due to the density threshold. The density profile at t = 3.30 is 
shown in Figure [7^. Again, at the center the density is lowered due to the continue nonlinear interaction of waves. 
Only ascendant features from last rebound are detectable due to the density threshold pn = 2.11. Figure [Tla and 
Figure [5^ have common sharp edges features indicating the presence of shock fronts. From Figure [T]; and Figure [5)3 
we see common regions formed by interaction of waves near the bases; the interaction of waves finish smoothing out 
the density profiles. The differences on the correspondent shapes are related with the differences in the triggering 
mechanisms (the spontaneous evolution vs. external forcing from the bases), however, the resulting cadences and 
characteristic speeds are the sa me as the observational ones. 

Zaitsev and Stepanov |l989l | introduced a model from where they derived relations between flare loop parameters 
as a function of loop length, period and amplitude of the oscillation: 

where T, n^, L, B are the temperature, density, length and magnetic strength of the loop. If we take r — 9.3 min, the 
dimensional characteristic time associated with the flare bursts; Q ~ 1, the quality factor or the number of individual 
flare bursts; and the intensity amplitude of the flare oscillating features as the numerical intensity amplitude, AH/ H, 
we can estimate the temperature, electron density and the magnetic field strength as 

T - 6.3 10^ K, Tie - 4.9 10^ cm'^, B - 8.7 G. 

To adjust the observations we performed several calculations varying the parameters AH/ < H >, Q and r respec- 
tively. We noted that the increasing of AH/ < H >, corresponds to more pronounced curves of the descending 
branches and that more complex patterns can be obtained depending on Q, the number of individual fiare burst cho- 
sen and the time distance between them, t. Taking into account the approximations we assumed and the numerical 
simplifications of the model we note that the values obtained are in good agreement with the observational values. 
The evaluation of the relative importance of the parameters and initial conditions to produce different patterns that 
could be compared to observations together with the consideration of nonidcal contributions will be a matter of a 
new work. 

A few words with respect to the consideration of diffusive terms in the equations must be added. Longitudinal 
modes can be affected by the thermodynamic properties of the medium. However, the fact that this type of modes 
is selectively seen in hot loops, where the mode damping is strongest, remains an open question. The scenario we 
are dealing with shows that nonlinear wave interaction processes are responsible of smoothing out sharp edges in a 
non-dissipative regime. This suggests that we are describing the stage where a non-dissipative transfer of energy from 
the large spatial scales to the small ones is predominantly at work before the small scales allow an efficient action of 
dissipation. 

A novelty of this work is that it proposes that the description in terms of shock waves is more common that 
what it was supposed. Low /3 value media have larger values of the Alfven speed than the sound speed. These two 
characteristic parameters of the coronal plasma determine the type of processes and energies associated to them that 
can be involved in its dynamics. The fact that in general, depending on the p value, the sound speed is the smallest 
one suggest that it must not be unusual that compressional perturbations reach shock wave speeds. Inhomogeneous 
features in the medium can easily couple the transverse alfvenic component with the longitudinal one and thus transfer 
energy to the acoustic mode allowing that hydrodynamic shock waves, guided by magnetic fields, can propagate. 
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FIG. 6: a) Non-dimensional descendant and ascendant density features as a function of the non-dimensional time, b) Non- 
dimensional descendant density features using the threshold H. x = 1 is equivalent to L = 150 Mm. One time division is 
equivalent to 18 min lOsec. 




FIG. 7: Non-dimensional density as function of the non-dimensional loop length taken from the model in subsection (4.2) at 
the a) non-dimensional time 2.23, two descending shock fronts 2:13 = 0.4; 2:23 = 0.6; ascending waves, corresponding to the 
rebound of the first shock fronts, with maxima at xi --^ 0.1; 12 ~ 0.9, both below the threshold ph ^ 2.11; b) non-dimensional 
time 2.57, two descending shock fronts at xis ~ 0.2; 2:23 ^ 0.8; ascending waves lower the density at the center of the figure; 
c) non-dimensional time 3.30, the interaction of waves smooth the profile, only the ascendant features are detected due to the 
density threshold, pn ~ 2.11. 



V. CONCLUSIONS 



We have integrated the MHD ideal equations of a slender flux tube to simulate the internal plasma dynamics of 
coronal post-flare loops. We could reproduced the observational behaviour of brightening features along magnetic 
threads of an event occurred on October 1st, 2001, i.e., characteristic speeds and times. Our intention was to analyze 
the plasma dynamics interior to postflare loops and its relation with flaring events, i.e. whether the oscillatory 
patterns could be associated to internal or external driving mechanisms. This is, whether non-dissipative second 
standing acoustic modes could be responsible for the induction of QPPs or if they are resulting features of an external 
forcing produced by the flaring event. We showed that, in the frame of our model, to reproduce the iteration of sliding 
down observations a time-dependent forcing from the basis is required and thus the interpretation that the QPPs are 
responses to the pulsational flaring activity is reinforced. 

We also found that high-speed downflow perturbations usually interpreted as slow magnetoacoustic waves can be 
better interpreted as slow magnetoacoustic shock waves. 



Case pi 


= p:i [cm' 


p2 (cm-^) 


Ti = n{K) T2{K) 


wave speed sound speed 


I 


10'' 


10^" 


10" 


10" 


0.81 


0.374 


II 


10^" 


10^" 


lO'' 


10" 


0.81 


0.374 


III 


10^^ 


10^^ 


10* 


lO'' 


0.81 


0.374 


IV 


10^^ 


10^^ 


lO'' 


10" 


0.81 


0.374 


V 


10^" 


10"' 


lO'' 


2.10" 


0.771 


0.265 


VI 


10^^ 


10^^ 


10"* 


2.10" 


0.771 


0.265 



TABLE I: Values of the six initial condition cases that adjust the observational data of the symmetric case (see Figure [T^). 
The wave front speed and the acoustic speed were nondimenzionalized with the Alfven speed {va = 1). 
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